NS Hemisphere Visualizations
KE Lotterhos
10/27/2020
- upload MPA locations
- link MPA locations to grid cells on globe
MPAs <- read.csv(paste0(data_dir, "notakeonly_coord.csv"))
head(MPAs)## Centroid_Longitude Centroid_Latitude
## 1 177.8167 -17.35000
## 2 -17.0500 20.78333
## 3 179.1040 -17.10960
## 4 179.9667 -16.85000
## 5 -159.7500 -18.86670
## 6 158.1257 6.79908tail(MPAs)## Centroid_Longitude Centroid_Latitude
## 304 8.559806 42.35926
## 305 3.164926 42.47173
## 306 6.390493 43.00589
## 307 -81.396041 19.39510
## 308 -111.091496 25.85519
## 309 -42.750000 -62.75000plot(MPAs$Centroid_Longitude, MPAs$Centroid_Latitude)
# Transform longitude from MPAs
MPAs$Centroid_Longitude_trans <- MPAs$Centroid_Longitude
MPAs$Centroid_Longitude_trans[which(MPAs$Centroid_Longitude<0)] <- 360 + MPAs$Centroid_Longitude[which(MPAs$Centroid_Longitude<0)]
dfresults <- read.csv(paste0(results_dir,"SigmaD_Hem.csv"))
dfresults <- dfresults[order(dfresults$No),]
head(dfresults)## No NN.sigma_2100_4.5_today NN.station_2100_4.5_today
## 17414 1 0.0008453082 30792
## 17415 2 0.0003516396 30792
## 17416 3 0.0018296424 30792
## 17417 4 0.0032007127 30792
## 17418 5 0.0074244659 30575
## 17419 6 0.0055364983 30575
## NN.Mdist_2100_4.5_today numPCs_2100_4.5_today NN.sigma_today_2100_4.5
## 17414 0.03673379 2 2.745675
## 17415 0.02368996 2 2.692024
## 17416 0.05405385 2 2.677252
## 17417 0.07151309 2 2.705799
## 17418 0.10900851 2 2.776901
## 17419 0.09409835 2 2.891746
## NN.station_today_2100_4.5 NN.Mdist_today_2100_4.5 numPCs_today_2100_4.5
## 17414 40066 3.196740 2
## 17415 40066 3.145595 2
## 17416 40066 3.131513 2
## 17417 40066 3.158725 2
## 17418 40066 3.226513 2
## 17419 40066 3.336049 2
## NN.sigma_2100_8.5_today NN.station_2100_8.5_today NN.Mdist_2100_8.5_today
## 17414 5.866236 30680 6.201382
## 17415 5.993259 30680 6.324892
## 17416 6.079694 30680 6.408981
## 17417 6.122091 30680 6.450239
## 17418 6.114727 30680 6.443073
## 17419 6.062258 30680 6.392015
## numPCs_2100_8.5_today NN.sigma_today_2100_8.5 NN.station_today_2100_8.5
## 17414 2 8.160708 40066
## 17415 2 8.292361 40066
## 17416 2 8.160708 40066
## 17417 2 8.209536 40066
## 17418 2 8.209536 40066
## 17419 2 8.292361 40066
## NN.Mdist_today_2100_8.5 numPCs_today_2100_8.5 lat long
## 17414 8.512412 2 -69.5 20.5
## 17415 8.470758 2 -68.5 20.5
## 17416 8.465380 2 -67.5 20.5
## 17417 8.501935 2 -66.5 20.5
## 17418 8.577215 2 -65.5 20.5
## 17419 8.693851 2 -64.5 20.5
## NN.sigma_1800_today NN.station_1800_today NN.Mdist_1800_today
## 17414 0.005340623 40066 0.0924152
## 17415 0.007501424 40009 0.1095737
## 17416 0.009629160 40067 0.1241975
## 17417 0.015490226 40009 0.1577081
## 17418 0.024201364 40066 0.1974682
## 17419 0.049860481 40066 0.2848815
## numPCs_1800_today NN.sigma_today_1800 NN.station_today_1800
## 17414 2 0.019575390 26505
## 17415 2 0.020981196 29981
## 17416 2 0.014089519 26082
## 17417 2 0.006465863 26082
## 17418 2 0.011645698 25749
## 17419 2 0.005112661 25973
## NN.Mdist_today_1800 numPCs_today_1800 lat_1800_today long_1800_today
## 17414 0.17743251 2 -69.5 379.5
## 17415 0.18374455 2 -69.5 378.5
## 17416 0.15036689 2 -68.5 379.5
## 17417 0.10170870 2 -69.5 378.5
## 17418 0.13663930 2 -69.5 379.5
## 17419 0.09041724 2 -69.5 379.5
## latshift_1800_today latshift lat_today_2100_8.5 long_today_2100_8.5
## 17414 NA 0 -69.5 379.5
## 17415 NA -1 -69.5 379.5
## 17416 NA -1 -69.5 379.5
## 17417 NA -3 -69.5 379.5
## 17418 NA -4 -69.5 379.5
## 17419 NA -5 -69.5 379.5
## lat_today_2100_4.5 long_today_2100_4.5
## 17414 -69.5 379.5
## 17415 -69.5 379.5
## 17416 -69.5 379.5
## 17417 -69.5 379.5
## 17418 -69.5 379.5
## 17419 -69.5 379.5tail(dfresults)## No NN.sigma_2100_4.5_today NN.station_2100_4.5_today
## 17408 40115 4.501010 554
## 17409 40116 4.525272 554
## 17410 40117 4.569456 554
## 17411 40118 4.606448 554
## 17412 40119 4.625980 554
## 17413 40120 4.632551 554
## NN.Mdist_2100_4.5_today numPCs_2100_4.5_today NN.sigma_today_2100_4.5
## 17408 4.879349 2 5.414412e-05
## 17409 4.902748 2 7.803501e-06
## 17410 4.945368 2 5.824597e-04
## 17411 4.981062 2 3.656808e-04
## 17412 4.999911 2 3.734047e-03
## 17413 5.006253 2 1.033971e-02
## NN.station_today_2100_4.5 NN.Mdist_today_2100_4.5 numPCs_today_2100_4.5
## 17408 8267 0.009295342 2
## 17409 24925 0.003528828 2
## 17410 12303 0.030490772 2
## 17411 21911 0.024158378 2
## 17412 16363 0.077249962 2
## 17413 16363 0.128716445 2
## NN.sigma_2100_8.5_today NN.station_2100_8.5_today NN.Mdist_2100_8.5_today
## 17408 8.292361 40112 9.377801
## 17409 8.292361 40113 9.465370
## 17410 8.292361 40113 9.475441
## 17411 8.292361 40113 9.455885
## 17412 8.292361 40113 9.425748
## 17413 8.292361 40113 9.404554
## numPCs_2100_8.5_today NN.sigma_today_2100_8.5 NN.station_today_2100_8.5
## 17408 2 0.9093877 9460
## 17409 2 1.2192752 8629
## 17410 2 1.4372294 8629
## 17411 2 1.6185109 8629
## 17412 2 1.7344557 8629
## 17413 2 1.7884077 8629
## NN.Mdist_today_2100_8.5 numPCs_today_2100_8.5 lat long
## 17408 1.423342 2 84.5 379.5
## 17409 1.733062 2 85.5 379.5
## 17410 1.945650 2 86.5 379.5
## 17411 2.120635 2 87.5 379.5
## 17412 2.231984 2 88.5 379.5
## 17413 2.283687 2 89.5 379.5
## NN.sigma_1800_today NN.station_1800_today NN.Mdist_1800_today
## 17408 0.0033782225 37019 0.07347198
## 17409 0.0041649559 37769 0.08159263
## 17410 0.0001037297 3229 0.01286606
## 17411 0.0047175381 7551 0.08684629
## 17412 0.0377153031 7551 0.24717272
## 17413 0.0646634655 7551 0.32537576
## numPCs_1800_today NN.sigma_today_1800 NN.station_today_1800
## 17408 2 1.728258 613
## 17409 2 1.931682 613
## 17410 2 2.170640 613
## 17411 2 2.350300 554
## 17412 2 2.424026 3223
## 17413 2 2.455741 3336
## NN.Mdist_today_1800 numPCs_today_1800 lat_1800_today long_1800_today
## 17408 2.226041 2 86.5 347.5
## 17409 2.420731 2 86.5 353.5
## 17410 2.648753 2 88.5 59.5
## 17411 2.819972 2 88.5 103.5
## 17412 2.890216 2 88.5 103.5
## 17413 2.920432 2 88.5 103.5
## latshift_1800_today latshift lat_today_2100_8.5 long_today_2100_8.5
## 17408 -2 -2 71.5 130.5
## 17409 -1 -1 73.5 119.5
## 17410 -2 -2 73.5 119.5
## 17411 -1 -1 73.5 119.5
## 17412 0 0 73.5 119.5
## 17413 1 1 73.5 119.5
## lat_today_2100_4.5 long_today_2100_4.5
## 17408 77.5 113.5
## 17409 83.5 239.5
## 17410 81.5 157.5
## 17411 79.5 218.5
## 17412 78.5 183.5
## 17413 78.5 183.5world <- map_data("world2")
ggplot()+
geom_polygon(data = world, aes(x=long, y = lat, group=group)) +
geom_point(data = MPAs, aes(x= (Centroid_Longitude_trans), y = Centroid_Latitude)) +
stat_summary_2d(data=dfresults, aes(x=long, y = lat, z= NN.sigma_today_2100_4.5), bins=80, alpha = 0.8)
# Transform longitude from this analysis
dfresults$long_trans <- dfresults$long
dfresults$long_trans[which(dfresults$long>360)] <- dfresults$long[which(dfresults$long>360)] - 360
# 380 should correspond to 20; 360 corresponds to 0
ggplot()+
# geom_polygon(data = world, aes(x=long, y = lat, group=group)) +
geom_point(data = MPAs, aes(x= (Centroid_Longitude_trans), y = Centroid_Latitude)) +
stat_summary_2d(data=dfresults, aes(x=long_trans, y = lat, z= NN.sigma_today_2100_4.5), bins=80, alpha = 0.8)
summary(MPAs$Centroid_Longitude) # - 180 to 180## Min. 1st Qu. Median Mean 3rd Qu. Max.
## -177.37 -121.95 -76.01 -10.41 143.40 179.97summary(dfresults$long) # 20 to 380## Min. 1st Qu. Median Mean 3rd Qu. Max.
## 20.5 138.5 207.5 208.8 295.5 379.5dfresults$MPA_nnIndex <- NA
dfresults$MPA <- FALSE
dfresults$MPA_nnDist <- NA
MPAs$No <- NA
for (i in 1:nrow(MPAs)){
q <- as.matrix(cbind(MPAs$Centroid_Latitude[i], MPAs$Centroid_Longitude_trans[i])) # lat, long
d <- as.matrix(dfresults[,c("lat", "long_trans")])
nn <- get.knnx(d, q, k=1)
MPAs$No[i] <- nn$nn.index
MPAs$nnDist[i] <- nn$nn.dist
MPAs$No_lat[i] <- dfresults$lat[nn$nn.index]
MPAs$No_long_trans[i] <- dfresults$long_trans[nn$nn.index]
#q
#dfresults[nn$nn.index,c("lat", "long")]
if(nn$nn.dist < 0.9){ # only include if station is close
dfresults$MPA_nnDist[nn$nn.index] <- nn$nn.dist
dfresults$MPA_nnIndex[nn$nn.index] <- nn$nn.index
dfresults$MPA[nn$nn.index] <- TRUE
}
}
# total number of MPAs included
length(which(MPAs$nnDist<0.9))## [1] 177# number of stations representing MPAs
sum(!is.na(dfresults$MPA_nnIndex))## [1] 94table(MPAs$No[which(MPAs$nnDist<0.9)])##
## 591 1098 1447 2657 4396 4856 8019 8091 8393 8452 8718 9503 9838
## 1 1 1 3 1 1 1 1 1 1 3 1 1
## 9839 10004 10088 10367 10369 10454 10530 10549 10550 10719 10825 10828 10830
## 3 1 1 1 1 1 3 1 3 2 1 1 2
## 10940 10946 11057 11059 11172 11180 11181 11188 11196 11439 11440 11441 11442
## 1 1 1 1 1 2 2 1 1 1 1 2 3
## 12168 12169 12272 12274 12534 12741 13594 13803 14377 15029 15157 15188 15189
## 5 2 1 2 2 1 1 1 1 3 1 1 1
## 15320 15350 15351 15885 16678 16833 17038 17197 18342 18765 18768 18920 18962
## 1 3 1 1 1 3 1 1 1 1 3 1 2
## 19095 19114 19266 19267 19571 23996 24133 24135 24263 24264 24390 24513 24634
## 2 5 2 2 1 1 1 1 6 10 7 4 3
## 24635 26464 27761 28671 28842 29090 29238 29512 29546 29665 29730 29731 29875
## 7 2 1 2 1 1 1 1 8 2 3 1 1
## 32070 33374 35810
## 1 1 1 #This should sum to 177
sum(table(MPAs$No[which(MPAs$nnDist<0.9)]))## [1] 177hist(dfresults$MPA_nnDist)
- look at degree of novelty for MPAs relative to non-MPA sites
ggplot(dfresults, aes(NN.sigma_today_2100_4.5,
fill = MPA)) +
geom_density(color="#e9ecef", alpha=0.3, position = 'identity') + theme_classic()
ggplot(dfresults, aes(NN.sigma_today_2100_8.5,
fill = MPA)) +
geom_density(color="#e9ecef", alpha=0.3, position = 'identity') + theme_classic()
- make arrow plots from 1800 to today to 2100 for MPAs
To Do 2: 1. Create fake environment envelope 2. Create fake locations with high ICV and different degrees of novelty 3. Create fake locations with low ICV and different degrees of novelty 4. Visualize
Visualize arrows
There is a confusing issue with this MPA analysis. I can’t just take the MPA column from dfresults and have it make sense all the time.
the 1800_today analysis takes a site today (given by lat,long) and compares it to the 1800 baseline. The NN in the 1800 baseline is lat_1800_today, long_1800_today. So if I want to draw an arrow from where a climate was in 1800 to where it is today today, I would draw from lat_1800_today to lat
the today_2100 analysis takes a site in 2100 (given by lat, long) and compares it to a today baseline. The NN in the today baseline is lat_today_2100, long_today_2100. So if I want to draw an arrow from where a climate is today in today’s baseline to where it will be in 2100, I would draw from lat_today_2100 to lat
ggplot() +
geom_polygon(data = world, aes(x=long, y = lat, group=group)) +
# Plot from 1800 to location today
# Arrows point from 1800 climate to MPA today
geom_segment(data = dfresults[dfresults$MPA,],
aes(x = long_1800_today,
y = lat_1800_today,
xend = long, yend = lat,
colour = NN.sigma_1800_today),
arrow = arrow(length = unit(0.1,"cm")))
ggplot() +
geom_polygon(data = world, aes(x=long, y = lat, group=group)) +
# Plot from today to 2100 4.5
# In 2100, what is the nearest neighbor today?
# Arrow points from today's environment (long_today_2100_4.5) to the MPA in the future (long)
geom_segment(data = dfresults[dfresults$MPA,],
aes(xend = long,
yend = lat,
x = long_today_2100_4.5,
y = lat_today_2100_4.5,
colour = NN.sigma_today_2100_4.5),
arrow = arrow(length = unit(0.2,"cm"))) 
ggplot() +
geom_polygon(data = world, aes(x=long, y = lat, group=group)) +
# Plot from today to 2100 8.5
geom_segment(data = dfresults[dfresults$MPA,],
aes(xend = long,
yend = lat,
x = long_today_2100_8.5,
y = lat_today_2100_8.5,
colour = NN.sigma_today_2100_8.5),
arrow = arrow(length = unit(0.2,"cm"))) 
Compare latitudes
poly <- data.frame(y=c(-90,90,90, -90, -90), x=c(-90,90,0,0,-90))
# If we consider a location in 2000 (lat), where is it's nearest neighbor 1800? (lat_1800_today)?
ggplot(dfresults) + geom_polygon(data=poly,aes(x,y), fill=adjustcolor("orange", 0.3)) +
geom_point(aes(y=lat, x=lat_1800_today,
color=NN.sigma_1800_today), alpha=0.2) +
scale_color_gradient2(low="red", high="blue", mid="grey", limits=c(0,8)) +
theme_classic() + ylab("Latitude in 2000") +
xlab("Latitude of nearest neighbor in 1800") + geom_abline(intercept=0,slope=1) 
# 2000 to 2100
# If we consider a location in 2100, where was the climate it came from in 2000?
ggplot(dfresults) + geom_polygon(data=poly,aes(x,y), fill=adjustcolor("orange", 0.3)) +
geom_point(aes(y=lat, x=lat_today_2100_4.5,
color=NN.sigma_today_2100_4.5),
alpha=0.2) +
scale_color_gradient2(low="red",
high="blue", mid="grey", limits=c(0,8)) +
theme_classic() + ylab("Latitude in 2100 RCP 4.5") +
xlab("Latitude of nearest neighbor in 2000") + geom_abline(intercept=0,slope=1) 
ggplot(dfresults) + geom_polygon(data=poly,aes(x,y), fill=adjustcolor("orange", 0.3)) +
geom_point(aes(y=lat, x=lat_today_2100_8.5,
color=NN.sigma_today_2100_8.5),
alpha=0.2) +
scale_color_gradient2(low="red",
high="blue", mid="grey", limits=c(0,8)) +
theme_classic() + ylab("Latitude in 2100 RCP 8.5") +
xlab("Latitude of nearest neighbor in 2000") + geom_abline(intercept=0,slope=1) 
Latitudinal shifts
## Future latitude - past latitude
## 1800 to today
## "lat" is latitude of query in 2000
## lat_1800 is latitude of NN in 1800
dfresults$latshift_1800_today <- NA
condN <- which(dfresults$lat>0)
dfresults$latshift_1800_today[condN] <- dfresults$lat[condN] - dfresults$lat_1800_today[condN]
# If 50 degrees in 1800 moved to 90 degrees in 2100 this
# would give 90-50 = 40
condS <- which(dfresults$lat<=0)
dfresults$latshift_1800_today[condS] <- -1*(dfresults$lat[condS]-dfresults$lat_1800_today[condS])
# If -50 degrees in 1800 moved to -90 degrees in 2100 this
# would give -1*(-90--50) = 40
## Today to 2100 RCP 4.5
## Future latitude - past latitude
dfresults$latshift_today_2100_4.5 <- NA
dfresults$latshift_today_2100_4.5[condN] <- dfresults$lat[condN] - dfresults$lat_today_2100_4.5[condN]
dfresults$latshift_today_2100_4.5[condS] <- -1*(dfresults$lat[condS] - dfresults$lat_today_2100_4.5[condS])
## Today to 2100 RCP 8.5
dfresults$latshift_today_2100_8.5 <- NA
dfresults$latshift_today_2100_8.5[condN] <- dfresults$lat[condN] - dfresults$lat_today_2100_8.5[condN]
dfresults$latshift_today_2100_8.5[condS] <- -1*( dfresults$lat[condS] - dfresults$lat_today_2100_8.5[condS])
par(mar=c(4,4,1,1), mfrow=c(3,1))
hist(dfresults$latshift_1800_today, col=adjustcolor("blue", 0.5), main="", xlab="Poleward shift from\nneareast neighbor in the past", breaks=seq(-90, 90, 1), xlim=c(-90,90))
hist(dfresults$latshift_today_2100_4.5,
col=adjustcolor("red", 0.2), main="", xlab="",
breaks=seq(-90, 90, 1), xlim=c(-90,90))
hist(dfresults$latshift_today_2100_8.5,
col=adjustcolor("green", 0.2), main="", xlab="",,
breaks=seq(-90, 90, 1), xlim=c(-90,90))
# color scale: if use "lat" it is the latitude in the later time
# color scale: if use "lat_1800_today" it is the latitude in the earlier time
# the x-axis is the shift of the latitude to that location
ggplot(dfresults) + geom_point(aes(x=latshift_1800_today, y = NN.sigma_1800_today, colour = lat), alpha=0.1)
ggplot(dfresults[dfresults$lat_1800_today>85,]) + geom_point(aes(x=latshift_1800_today, y = NN.sigma_1800_today, colour = lat), alpha=0.1)
# high latitudes in 1800 can only stay in the same place or move southward - color is latitude in 2000
ggplot(dfresults[dfresults$lat_1800_today>70 & dfresults$lat_1800_today<85,]) + geom_point(aes(x=latshift_1800_today, y = NN.sigma_1800_today, colour = lat), alpha=0.1)
ggplot(dfresults[dfresults$lat_1800_today<60 & dfresults$lat_1800_today>40,]) + geom_point(aes(x=latshift_1800_today, y = NN.sigma_1800_today, colour = lat_1800_today), alpha=0.1)
# 40-60 degress latitude in 1800 more northward
ggplot(dfresults) + geom_point(aes(x=latshift_today_2100_4.5, y = NN.sigma_today_2100_4.5, colour = lat_today_2100_4.5), alpha=0.1)
# colored by latitude in 2000, shows which way it shifted
ggplot(dfresults) + geom_point(aes(x=latshift_today_2100_8.5, y = NN.sigma_today_2100_8.5, colour = lat_today_2100_4.5), alpha=0.1)
# The RCP 8.5 shows that there aren't a lot of poleward shifts - because climates are so novel, the NN calculation becomes meaningless
# I wonder if some of the really negative latitudinal shifts are driven by distortion of the PC axes caused by the data not being multivariate normal.